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ABSTRACT 

The lowest order sigma-transformed momentum equation given by Mellor (J. Phys. Oceangr. 2003) 
takes into account a phase-averaged wave forcing based on Airy wave theory. This equation is shown to 
be generally inconsistent due to inadequate approximations of the wave motion. Indeed the evaluation 
of the vertical flux of momentum requires an estimation of the pressure p and coordinate transformation 
function s to first order in parameters that define the large scale evolution of the wave field, such as the 
bottom slope. Unfortunately there is no analytical expression for p and s at that order. A numerical 
correction method is thus proposed and verified. Alternative coordinate transforms that allow a sepa- 
ration of wave and mean flow momenta do not suffer from this inconsistency nor require a numerical 
estimation of the wave forcing. Indeed, the problematic vertical flux is part of the wave momentum flux, 
thus distinct from the mean flow momentum flux, and not directly relevant to the mean flow evolution. 



1. Introduction 

D Wave-induced motions are of prime importance in 
the upper ocean, and in the coastal ocean (e.g. Ard- 
huin et al. 2005 for a recent review). Therefore, usual 
three-dimensional primitive equations must be modified 
to account for waves. Among such modified equations, 
those based on surface-following coordinates provide 
physically sound definitions of velocities right up to the 
free surface, allowing a proper representation of surface 
shears and mixing on a vertical scale smaller than the 
wave height (i.e. a few meters). Any change of coor- 
dinate adds some complexity in the derivation, but the 
final equations can be relatively simple because part of 
the advective fluxes are removed, and boundary condi- 
tions may be simplified. A new set of such equations 
was recently derived by Mellor (2003) using a change of 
the vertical coordinate only, arguably the simplest pos- 
sible. Mellor's (2003) set of equations was originally 
derived for monochromatic waves, but it is easily ex- 
tended to random waves (e.g. Ardhuin et al. 2004, eq. 
8). Unfortunately, we show here that these equations, 
in the form given by Mellor, are not consistent in the 
simple case of shoaling waves without energy dissipa- 
tion. A modification is proposed to solve the problem, 
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but it requires a numerical evaluation of the wave forc- 
ing terms. This difficulty is due to the choice of aver- 
aging, and the same problem arises with the alternative 
Generalized Lagrangian Mean equations of Andrews and 
Mclntyre (1978, eq. 8.7a, hereinafter aGLM). Both Mel- 
lor's and aGLM equations describe the evolution of a 
momentum quantity that contains the three-dimensional 
wave (pseudo)-momentum (hereinafter called 'wave mo- 
mentum' for simplicity, see Mclntyre 1981 for details). 
Writing an evolution equation for this quantity requires 
an explicit description of the complex vertical fluxes of 
wave momentum that are necessary to maintain the ver- 
tical structure of the wave field in the surface gravity 
waveguide. 

2. The problem: wave motions and wave-following 
vertical coordinates 

We discuss here the simple case of monochromatic 
waves of amplitude a and wavenumber k propagating in 
the horizontal x direction, with all quantities uniform in 
the other horizontal y direction. The surface and bottom 
elevations are r)(x, t) and —h(x), respectively, so that the 

local mean water depth is D(x, t) = h(x) + rj{x, t), with 
the overbar denoting an Eulerian average over the wave 
phase. We shall assume that the maximum surface slope 
is a small parameter E\ = ka <C 1, and that the Eule- 
rian mean current u in the x-direction is uniform over the 
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depth. Thus uj will denote the radian wave frequency re- 
lated to k by the linear wave dispersion relation (e.g. Mei 
1989), 



ku 



ku+ [gk tanh(fcZ))] 



1/2 



(1) 



Finally, we assume that the water depth, current and 
wave amplitude change slowly along the x-axis with 
a slowness measured by a second small parameter e 2 
taken to be the maximum bottom slope. We thus assume 

\{dD/dx)\ < e 2 , \{da/dx)\ < e 2 , \{du/dx)/{a)\ < 



e 2 , \{da/dt)/{aa)\ 



< 



e 2 , \k(du/dt)/a 2 



< 



?2, 



(dD /dt) k/<j\ < £2 . The conditions on the bottom slope 
and current gradients are consistent with the condition on 
the wave amplitude gradient because in steady conditions 
the wave amplitude would change due to shoaling over 
the current and/or bottom. 

The vertical coordinate z is implicitly transformed into 
Mellor's c; coordinate through 



z = s (x, <;, t) = T) + <;D + '. 



(2) 



with J defined by Mellor's eq. (23b) as 

s = So = aFss cos (kx — Lut) (3) 

and the vertical profile function F$ s defined by 

_ sinh [kD(l + Q] _ sinh [k(z + h)} a_ 
55 " sinh (kD) ~~ sinh (kD) + ( D>' 

(4) 

The coordinate transformation from z to c; has the very 
nice property of following the vertical wave-induced mo- 
tion, at least for linear waves on a flat bottom, and to 
first order in e\. In that case the iso-c; surfaces are ma- 
terial surfaces, and the fluxes of horizontal momentum 
through one of these surfaces are simply correlations of 
pressure p times the slope of that surface ds jdx (figure 
1 .c), which replaces the wave-induced advective flux tiro 
in an Eulerian point of view (figure l.a). More generally, 
when averaging is performed following water particles 
over their trajectory (Lagrangian) or over their vertical 
displacement (Mellor-sigma), the corresponding advec- 
tive flux of momentum UiUj is replaced by a modified 
pressure force (figure 1). Q 

Using his coordinate transform, Mellor (2003) ob- 
tained a phase-averaged equation for the drift current 
U = u + us where us is the Stokes drift, i.e. the mean 
velocity of water particles induced by fast wave-induced 
motions. U is strictly defined as the phase-average parti- 
cle drift velocity when following the up-and-down wave 
motion, and u = U — us is a quasi-Eulerian mean current 
(Jenkins 1986, 1987). Below the wave crests u is equal, 
to second order in the wave slope, with the Eulerian mean 
current u (figure 2). 



'For the Generalized Lagrangian Mean (GLM) only the contribu- 
tions to lowest order in s\ are indicated. Indeed, in GLM the wave- 
induced advective flux is not strictly zero, but of higher order, since 
the average only follows a zero-mean displacement with a residual ad- 
vection, contrary to a truly Lagrangian mean with zero advection (e.g. 
Jenkins 1986). 
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FIG. 1 . Wave-induced fluxes of horizontal momentum in Eulerian, Gen- 
eralized Lagrangian, and Mellor-sigma averages of the flow. Viscous or 
turbulent fluxes are neglected for simplicity. Distorted squares repre- 
sent an elementary fluid volume and its position at four phases of the 
wave cycle, and the large arrow indicate the local wave orbital veloc- 
ity. The horizontal and vertical fluxes of the horizontal momentum are 
represented by smaller arrows. Their expression are given to lowest or- 
der, without Jacobian corrections due to a change of volume (e.g. this 
results in u 2 + p becoming Mellor's u 2 + pds/ds). 



Mellor's horizontal mean momentum equation (34a) is 
reproduced here for completeness, in our conditions with 
a flow restricted to the vertical x, z plane, a constant water 
density, no Coriolis force, and no turbulent fluxes and the 
atmospheric mean pressure set to zero (wind-wave gener- 
ation due to air pressure fluctuations is absorbed in F x s), 



dDU dDU 2 dflU 



dt 



dx 



gD 



drj 
dx 



F xx + F x3 . (5) 
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FIG. 2. Averaging procedures (left) and examples of resulting velocity 
profiles (right) in the case of (a) Eulerian averages, (b) the Generalized 
Lagrangian Mean, (c) sigma transforms (Mellor 2003). The thick black 
bars connect the fixed points (;r, z) where the average field is evaluated, 
to the displaced points (x, z) + (£1, £3) where the instantaneous field 
is evaluated. For averages in moving coordinates the points (x, z) + 
(^1,^3) at a given vertical level £ are along the gray lines. The drift 
velocity is the sum of the (quasi-Eulerian) current and the wave-induced 
mass transport. In the present illustration an Airy wave of amplitude 3 m 
and wavelength 100 m in 30 m depth, is superimposed on a hypothetical 
current of velocity u(z) = —0.5 — 0.01 z m/s for all 2 < (,(x). The 
quasi-Eulerian current profile is not represented in (c) since it is not 
directly given in Mellor's theory, although it can obviously be obtained 
by taking the difference of the other two profiles. 



On the right-hand the first other term 

dS xx 



Fr 



dx 



d I „„ _<9s 
— Du 2 + P — 

OX \ OS 



(6) 



represents the convergence of a horizontal flux of hori- 
zontal momentum that accelerates the mean drift velocity 
U. 

The other term 



as 



x 3 



a:3 



— (pds/dx 



(7) 



represents a similar convergence of a vertical flux of hor- 
izontal momentum. 

Defining g as the acceleration due to the apparent grav- 
ity, p and s are of the order of ga and a respectively. In 
general p and s are almost in phase, thus the flux S X 3 is 



of the order of ga da/dx, and the force F X 3 is of the or- 
der of gDka da/dx. Thus, in the case of shoaling waves, 
F X 3 is of the order of gDe\e2. 

Mellor estimated the vertical momentum flux S X 3 from 
(fj]) and the corresponding lowest order wave-induced 
kinematic pressure on ? levelfl 

p = Po = ga(F C c - F S s)cos(kx - u>t), (8) 

g is the acceleration due to the apparent gravity, and the 
vertical profile function Fee is defined by 



F, 



cc 



cosh [£;£>(! +?)] cosh [k(z + h)] 



cosh (kD) 



cosh (kD) 



-0( 



D 
(9) 

For non-dissipating shoaling waves, the right hand side 
terms of eq. <(5j are of order gDe\ei- The estima- 
tion of ^3 thus requires the knowledge of pj (gD) and 
ks to order £162, for which Airy theory is insufficient. 
In particular, this estimation demands a formal defini- 
tion of s", not given by Mellor (2003). Further, eq. (7) 
is only valid if the wave-induced velocity fl through q 

levels is zero, or at least, yields a negligible flux Slu 
and a negligible mean Jacobian-weighted vertical veloc- 
ity £1 = £1/(1 + ds/ds/D). This is not the case over a 
sloping bottom with Mellor's (2003) s function. 

a. Formal definition of the the coordinate change s 

For a general surface <; defined implicitly by z = 
s(x, <;,£), the ? velocity component £1 is (e.g. Mellor 
2003 eq. 20), 



£2 



d(z — s) d(z — s) 



dt 



dx 



dt 
dx 



ds_ 

dt' 



(10) 



with s = rj + <;D. 

In the spirit of Mellor's (2003) derivation, the <; lev- 
els should be material surfaces for wave-only motions, 
so that one may neglect the vertical flux of momentum 

(U + u) n/(l + ds/ds/D). 

Using the wave-induced vertical and horizontal 
displacements, £3(3;,^, t) and £i(x, S,t), defined by 
d^i/dt = Ui(x\ + £i,z + £3,*), we redefine the wave 
part of s, 

d~s 



= £3 



1 dx 



(ID 



The first term £3 corresponds to Mellor's definition while 
the second is a O(ea) relative correction. This definition 
yields a wave-induced vertical velocity Q = —udl'/dx 
through the iso-<j surfaces redefined by z = s (x, t) = 
rj + qD + s'. If u <C u, as in the examples below, then 



2 This pressure includes a hydrostatic correction due to the vertical 
displacement. 



4 



JOURNAL OF PHYSICAL OCEANOGRAPHY 



Volume 



is of a higher order compared to that given by Mellor's 
(2003)s(eq.|2j. 

b. Wave-induced vertical displacements and pressure 
over a sloping bottom 

A WKBJ approximation using Airy's theory is suffi- 
cient for estimating dF xx /dx because the horizontal gra- 
dient of any wave-averaged quantity </> is of order £2^0. 
On the contrary, the other force F X 3 is affected by mod- 
ifications s[ and pi to the local-fiat-bottom solutions s', 
and p. 

For small bottom slopes, s[ and pi are expected to be 
of the order of E2S 1 and E2P, i.e. of order 0,62 and gae\£i, 
respectively. Thus ffs±/dx is of order \1ae2 — £i£2, and 
is expected to be in phase with the wave-induced pressure 
dHJ, of order ga, giving another term of order gDe\e2 
omitted by Mellor in his estimation of d(pds/dx)/d<,. 
The modification of the pressure can be obtained from 
the modification of the velocity potential, and it may be 
in phase with 8Jq /dx, thus also contributing at the same 
order to F x %. 

In order to be convinced of the problem, one may con- 
sider the case of steady monochromatic shoaling waves 
over a slope without bottom friction, viscosity or any kind 
of surface stress. We also neglect the Coriolis force. In 
this mathematical experiment, the flow is purely irrota- 
tional. We consider that the non-dimensional depth kH 
is of order 1, and that there is no net mass flux across 
any vertical section. In that case the mean current and the 
Stokes drift are of the same order, i.e. of the order Ce\ 
with C the phase speed. The mean current exactly com- 
pensates the divergence of the wave-induced mass trans- 
port, and the mean sea level is lower in the area where the 
wave height is increased (Longuet-Higgins 1967) 



rj(x) = - 



kE 



k Eo 



sinh(2fcD) sinh(2fc -D ) 



(12) 



where the subscript correspond to quantities evaluated 
at the offshore boundary of the domain. 

Since wave forcing is steady, the Eulerian mean cur- 
rent response is steady (e.g. Rivero and Sanchez- Arcilla 
1994, McWilliams et al. 2004, Lane et al. 2006), and 
thus the Lagrangian mean current is also steady. Thus 
the first term in (0 is zero and the second is of order 
DC 2 e\ 62/ D ~ gDe\e2- The vertical mean velocity 51 
can be estimated from the steady mass conservation equa- 
tion, 

dDU dn n 

-t— + ^t = Q (13) 

ox <X 

where the first term is of order DCe\e2j 'D and the sec- 
ond is of order ft. Thus the third term in <(5j is of order 
C 2 e\e2 — gDe\e2- The remaining terms in (0 are of 
order £^£2, giving the lowest order momentum balance 



Feta Fxx Fx3 

d dS d 

~ D -£- (9V) k~ + Trpds/dx= 0. 

ox Ox o<; 



For reference the corresponding lowest order Eulerian 
mean balance is (e.g. Rivero and Sanchez-Arcilla 1994, 
Lane et al. 2006) 



d ( - — 



du 2 
dx 



duw 
~dz~ 



= 0, 



(15) 



Only the hydrostatic pressure gradient is present in both 
the Eulerian and Mellor-sigma balances, because the 
other terms represent a different balance, including wave 
momentum in the latter (see figure 2). 

Equation ( fT4l is now tested numerically. We take a 
Roseau-type bottom profile (1976) defined by x and z 
coordinates given by the real and imaginary part of the 
complex function 



Z(x') 



hx(x' - ia) + (h 2 - hi) ln(l + e x ~ iQ ) 



(14) 



(16) 

With a = 15vr/180, hi = 6 m and h 2 = 4 m (fig- 
ure 1), and a radian frequency uj = 1.2 rad s _1 (i.e. a 
frequency / = 0.2 Hz), the non-dimensional water depth 
varies between 0.85 < kH < 1.1. The reflection co- 
efficient for the wave amplitude is 1.4 x 10 -9 (Roseau 
1976), so that reflected waves may be neglected in the 
momentum balance. We illustrate the force balance ob- 
tained for waves with an offshore amplitude ao = 0.12 m, 
which corresponds to a maximum steepness £ X = ka — 
2.6 x 10 -2 equal to the maximum bottom slope £ 2 = £1. 
The change in wave amplitude is given by the conser- 
vation of the wave energy flux (see Ardhuin 2006 for 
a thorough discussion), and the wave phase ip is taken 
as the integral over x of the local wavenumber, so that 
dip/dx — k. The various terms are then estimated us- 
ing second order finite differences on a regular grid in <; 
coordinates, with 201 by 401 points covering the domain 
shown in figure 3. a. The three terms in eq. ( fl4] > are shown 
in figure 3. 

We have verified that the depth-integrated forces are in 
balance, within 0.1% of F cta . However, at most water 
depths there is a large imbalance, of the order of the in- 
dividual forces ( i.e. gDefe 2 ) ,up to 180% of F eta . This 
contradicts the known steady balance obtained from the 
Eulerian-mean analysis of Rivero and Sanchez-Arcilla 
(1994). 

For the case of shoaling waves without breaking the 
three-dimensional equations of motion of Mellor (2003) 
are not consistent to their dominant order, because of an 
improper approximation of S X 3 . This conclusion holds 
for any relative magnitude of the wave and bottom slopes 
£1 and £2. 

c. Wind-forced waves 

Clearly, any deviation of the wave-induced fields s', 
p, and u from Airy-wave theory may have strong effects 
on the vertical momentum flux term S x s. Another ex- 
ample of such a situation, correctly described by Mellor, 
is the case of wind-wave generation. We briefly address 
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FIG. 3. (a) Snapshot of the pressure field for a slo wly varying Airy 
wave ov er a the bottom topography given by eq. (16). The forces in the 
balance 1141 are shown in panels b, c and d, with their sum in panel e, 
all estimated from Mellor's analytical expressions. All forces have been 
normalized by gDe s . N.B. in the case shown here ei = £2 = £■ 



it here because the full solution can has not been given 
previously. Mellor focused on the wind-wave generation 
contribution to the vertical momentum flux term pds/ dx 
term. This equals the wave-supported wind stress at the 
sea surface, and, below, it explains the growth of the wave 
momentum profile with the same profile as that of the 
Stokes drift (Mellor 2003). 

In horizontally uniform conditions, the wave amplitude 
is a function of time only, and for the sake of simplicity 
we shall solve the problem in the frame of reference mov- 
ing at the velocity at which the wave phase is advected by 
the current. We write the wave-induced non-hydrostatic 
kinematic Eulerian pressure in the formp^ = Peo+Pew, 
the elevation as 77 = £ + ( w and the velocity potential 
as = (f> + 4> w , in which the subscript refers to the 
primary waves, and the w subscript refers to the added 
components in the presence of wind forcing. Taking a 
primary surface elevation of the form (q — a cos tp with 
the the phase -0 = kx — at, Mellor considered an atmo- 
spheric kinematic pressure fluctuation in quadrature with 
the primary waves 



Pa 



-gf3—as'mifj, 

Pa 



(17) 



with /3 a small non-dimensional wave growth factor, and 
p w and p a the densities of water and air respectively. He 
then assumed that the water-side wave-induced pressure 
was of the form 



PMellor 



-g(3a 



cosh [k(z + h )] 
cosh (kD) 



sin-0. (18) 



Implicitly V w is zero, and for his purpose <j> w was irrel- 
evant. We shall now also determine S w . The continuity 
of dynamic pressures at the surface i^J 



PEv 



■9S V 



- g (3a sin ip at z — Q. 



(19) 



A solution is obtained by solving Laplace's equation 
with proper boundary conditions, to first order in (3. The 
boundary conditions include the Bernoulli equation, 



9<j> . Pa 
"57 = ~9C Pa, 

at p w 



at 



(20) 



in which non-linear terms have been neglected because 
they are the sum of products of the form V0o ■ V£o, un- 
changed from the case without wind, and terms of the 
form V<f) w • V£o, which are negligible compared to the 
left-hand side terms for primary waves of small slope. 
Similarly, the surface kinematic boundary condition is 
linearized as 



at 



dz dt 
The combination of both yields 

dt 2 



C- 



8(j) p a dp a - 

= kt at z = C- 

oz p w at 



(21) 



(22) 



4> w is also a solution of Laplace's equation with the bot- 
tom boundary condition dcf) w /dz — at z = —h. With 
the fully resonant atmospheric pressure ( fTTI i envisaged by 
Mellor, one has 



Co 


= a(t)costjj, 


(23) 


00 


= F cc smV, 
a 


(24) 


Peo 


= ga(t)F cc cos ip, 


(25) 


da(t) 
dt 


/3aa(t) 
2 ' 


(26) 


4> w 


= f3g—F cc [Acosip- 


hflcos^'] (27) 



with ip' = (kx + at). The elevation and under-water 
non-hydrostatic pressure corresponding to cf> w are given 
by (|2"TT i and the linearized Bernoulli equation 



8_ 

at 



-Peo - PEv 



(28) 



3 Here the pressure is Eulerian. For correspondance to Mellor's pres- 
sures on f levels, one should take p = pe — gs 1 ■ 
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yielding 



(1999) have used 



Cte = £-[(l-,4)sinV> + Bsin^'] (29) 

PEw = gP-^Fcc [— (1 + A) sin^ + Bsmij)'] . 

(30) 

Mellor's expression for p w (eq. 18) is obtained by re- 
placing ( w and pew in ( fT9l , giving A = 1. One may 
take S = 1 to have lu = at i = 0, or more simply 
B = 0, which gives ( w = 0, and p_B„, = F C cPaPa/Pw 
The choice of B has no dynamical effect. In the present 
case 4> w should give a contribution to S xx because it is 
in phase with </>o, but this is a relative correction of order 
/3, thus negligible. To the contrary, the contribution of 
p w to (pds/dx) is quite important, because for uniform 
horizontal conditions this flux is otherwise zero. 

3. A solution to the problem ? 

Contrary to that particular wind-forcing term, there is 
no simple asymptotically analytical correction for p and 
s' that can account for the bottom slope and wave field 
gradient. A major problem in this situation is that the 
wave velocity potential becomes a non-local function of 
the water depth. The velocity potential and pressure fields 
may only be investigated analytically for plane beds (e.g. 
Ehrenmark 2005) or specific bottom profiles. Numeri- 
cal solutions for the three-dimensional wave motion are 
generally found as infinite series of modes (e.g. Mas- 
sel 1993). The velocity potential for any of these modes 
satisfies Laplace's equation with a local vertical profile 
F„ proportional to cos(fc„z + k n h) and a dispersion re- 
lation a 2 = gk n t&n(k n D). The local amplitudes of 
these modes are non-local functions of the water depth, 
and may be obtained numerically with a coupled-mode 
model (Massel 1993). This non-local dependance of the 
wave amplitude on the water depth arises from the elliptic 
nature of Laplace's equation, satisfied by the velocity po- 
tential in irrotational conditions. The series of modes can 
be made to converge faster by adding a 'sloping bottom 
mode' that often accounts for a large part of the correc- 
tion and is a local function of the depth and bottom slope. 
It is thus of interest to see if that correction only, with- 
out the infinite series, may provide a first order analytical 
correction to Mellor's momentum flux S X 3. 

Following Athanassoulis and Belibassakis (1999), one 
may define the velocity potential for that mode as 



dh 
dx 



aaDF(z) cos ip, 



(31) 



In order to satisfy the bottom boundary condition w — 
d(pi/dz = —^d<fio/dx, the function F should ver- 
ify DdF(— h)/dz = 1/ sinh(kH) and the satisfaction 
of the surface boundary condition may be obtained with 
F(0) = dF(0)/dz = 0. Athanassoulis and Belibassakis 



F = F AB 



smh(kD) 



3 

D 



D 



and Chandrasekera and Cheung (2001) have used 



(32) 



F = F, 



1 



Ch 



fcDsinli (kD) 



[l — cosh(fcz — fc£)] 



(33) 

With these choices cf>i does not satisfy exactly 
Laplace's equation, and thus requires further corrections 
in the form of evanescent modes. An infinite number 
of other choices is available, either satisfying Laplace's 
equation or the surface boundary conditions, but never 
both, so that each of these solutions is only approxi- 
mate, and the exact solution is given by the infinite series 
of modes, which can be computed numerically for any 
bottom topography (e.g. Athanassoulis and Belibassakis 
1999, Belibassakis et al. 2001, Magne et al. 2006). 

The vertical displacement and Eulerian pressure cor- 
rections are given by time integration of the vertical ve- 
locity and the linearized Bernoulli equation, 



t dh n dF • / 

ax az 
dh 

Pei — —aDF(z) sin w 
dx 



(34) 
(35) 



Thus, in absence of wind forcing but taking into ac- 
count the 'sloping bottom mode' to first order in the bot- 
tom slope, the wave-induced flux of momentum through 
iso-? surfaces is 



P 



JS3_ 
dx 



gka 2 dh 
2 dx 



(Fee - 

qka 2 dh 
+ 

2 da; 

-FscF + Fss 



F 



ss) 



ga d(aF S s) 



D- 



dF 



D- 



dz 
dF 
d^ " 



dx 
h<;F cs 



(36) 



with F C s = cosh [kD( l+<; )]/sinh (kD). The first 
line is the term given by Mellor (2003). The second 
line arises from the correction due to the difference be- 
tween s' and s, and the third line arises due to corrections 
Pi = Pei — gs'i to the pressure on ^ levels. These ad- 
ditional term are of the same order as the first term, and 
have no flux at the bottom and surface. Thus the depth- 
integrated equations including that term also comply with 
known depth-integrated equations (e.g. Smith 2006). 

In the case chosen here Fch gives a net momentum bal- 
ance closer to zero than Mellor's (2003) original expres- 
sion (figure 4). However, the remaining error is signifi- 
cant. Thus one cannot use only that mode, and the con- 
tribution of the evanescent modes have to be computed, 
which can only be done numerically. 
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J(Feta + Fxx + Fx3) dx /jFeta dx 



FIG. 4. Net forces in the momentum balance 1141 for steady shoaling 
waves over a smooth bottom profile. The net force has been integrated 
over x and normalized by a similar integration of the the hydrostatic 
pressure force F e t a . Several solutions are obtained. One corresponds 
to Mellor's orginal expression, one possible analytical correction using 
i<Ch> an d numerical estimations using the NTUA model, with various 
numbers of modes. 



A numerical evaluation of the forces was performed 
using the NTUA model (Athanassoulis and Belibassakis 
1999). The NTUA solution was obtained in a domain 
with 401 points in the horizontal dimension. For the 
small bottom slope used here, the model contains a nu- 
merical reflection R = 0.002 much larger than the analyt- 
ical value given by Roseau (1976). However, this only in- 
troduces a modulation, in the x direction, of the estimated 
forces. This modulation is significant but still relatively 
smaller than the average. The net force estimated from 
NTUA results is found to converge to the expected force 
balance described by eq. (fl4l as the number of evanes- 
cent modes is increased (figure 4). In this calculation the 
values of F xx do not differ significantly from those esti- 
mated using Mellor's analytical expressions, as expected. 
The only significant difference between the NTUA nu- 
merical result with 10 modes, and Mellor's analytical ex- 
pression is found in F X 3, with a much stronger value near 
the surface in the numerical result, allowing a balance 
with the strongly sheared F xx (figure 4). 

4. Conclusions 

Mellor (2003) changed the vertical coordinate from 
C to z, using an implicit function s in two parts, z = 
s(x, y,(,,t) + s(x,y,<;,t) with s changing only slowly 
in space and time and s representing the faster wave- 
induced change of vertical coordinate. If the ? levels 
are material surfaces, then the momentum flux S X 3 = 
p^ds/dx is the surface-following coordinate counterpart 
of the Eulerian vertical momentum flux term uw dis- 
cussed by Rivero and Arcilla (1995), with p£ the wave- 



induced pressure at the displaced position (in the surface- 
following coordinates). However, pfids 1 /dx and uw do 
not represent the same physical quantity since the former 
contains wave momentum, which is not included in the 
latter. 

Just like the Eulerian momentum flux uw is modified 
by the bottom slope, wave amplitude gradients, wind- 
wave generation, boundary layers, or vertical current 
shears, these effects also modify S x a. But in these situa- 
tions, the <j levels as defined by Mellor (2003) are not ma- 
terial surfaces, and a missing Eulerian-like flux term ilw 
would have to be added to correct the momentum equa- 
tions, with fl the wave-induced velocity across s levels. 
Alternatively, we propose to replace s with s', defined by 
eq. ( fTTT ) such that <r levels are closer to material surfaces, 
i.e. so that fl is of a higher order. 

Whether the original J or our corrected s 7 is used, the 
wave-induced momentum flux S X 3 must be estimated to 
first order in the bottom slope £2 for consistency. This 
requires an Oiez) estimation of both p$ and J or 1'. Un- 
fortunately there is no analytical 0(^2) expression for the 
wave motion. Thus Mellor's equations, even when cor- 
rected, require a computer-intensive solution that is gen- 
erally not feasible. For example, Magne et al. (2006) 
only included a total of five modes in their calculation of 
wave propagation over a submarine canyon. In an exam- 
ple shown here, this small number of modes is insufficient 
for an accurate estimation of wave-forcing terms. 

The trouble with these equations can be avoided by us- 
ing, instead, equations of motion for the quasi-Eulerian 
velocity u = U - u s (Jenkins 1986, 1987, 1989). Such 
equations have been obtained in the limit of vanishing 
wave amplitude using an analytical continuation (e.g. us- 
ing a Taylor expansion) of the current profile across the 
surface (Mc Williams et al. 2004). A general and explicit 
solution can also be obtained from the exact General- 
ized Lagrangian Mean (GLM) equations of Andrews and 
Mclntyre (1978a) expanded to second order in the surface 
slope £1 (Ardhuin et al., manuscript submitted to Ocean 
Modelling, see |http://arxiv.org/abs/physics/0702067|). 
In these, the equation for the horizontal quasi-Eulerian 
momentum involves no flux term like pds/dx because 
this corresponds to the flux p^d^/ dx\l + 0(e)] of wave 
momentum us (Andrews and Mclntyre 1978b, eq. 2.7b), 
not directly relevant to the problem of mean flow evolu- 
tion (see also Jenkins and Ardhuin 2004). This flux of 
wave momentum only appears in evolution equations for 
the total momentum U, such as given by Mellor (2003), 
or the 'alternative' form of the GLM equations (Andrews 
and Mclntyre 1978, eq. 8.7a). 

For that reason, the equations for the quasi-Eulerian 
velocity u are simple and consistent in their adiabatic 
form (without wave dissipation), at least to lowest order 
in wave slope and current vertical shear, for which ana- 
lytical expressions exist for the wave forcing terms. Fur- 
ther details on the relationships between all these equa- 
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tions, and further validation against numerical solutions 
of Laplace's equation can be found in Ardhuin et al. (sub- 
mitted manuscript). 
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